start;
afa = 9.*(10.^(-3:0.5:1))*glv.deg;  % half-apex angle
f = 4;  w = 2*pi*f;                 % frequency

res1s = zeros(8,13,5);
for ss = 2:6
[nn, ts, nts] = nnts(ss, 0.01);  % subsamples & sampling interval
T = 1;                           % simulation time
len = T/ts; res = zeros(fix(len/nn), 3*4+1);
resCA = zeros(fix(len/nn), size(res,2), 10);

for ii = 1:length(afa)
[wm, qr] = conesimu(afa(ii), f, ts, T);
q1 = qr(1,:)'; q2 = q1; q3 = q1; q4 = q1; q5 = q1;  ki = 1;
for k=1:nn:len-nn+1
    k1 = k+nn-1;
	wmi = wm(k:k1, :);
	q0 = qr(k1+1,:)';
	phim1 = cnscl(wmi, 0);      q1 = qmul(q1,rv2q(phim1));  % optimal method
    phim2 = cnscl(wmi, 1);      q2 = qmul(q2,rv2q(phim2));  % polynomial method
    dQ = GetdQ_CP(wmi,ts,ss*2); q3 = qmul(q3,dQ);           % LP
	res(ki,:) = [q2att(q0); qq2phi(q1,q0); qq2phi(q2,q0); qq2phi(q3,q0); k*ts]';  ki = ki+1;
end
t = (1:size(res,1))'*nts;
[ebx, dphim] = conedrift(nn, ts, afa(ii), f);       % coning drift error

res1s(ii,:,ss-1) = [abs(res(end,1:12))./glv.dph afa(ii)./glv.deg];
resCA(:,:,ii) = res;
end
end

myfigure('1',0);
loglog(res1s(:,end,1), res1s(:, 4,1), 'b--', 'Marker','^', 'linewidth',1.5); grid on; hold on;
loglog(res1s(:,end,1), res1s(:, 7,1), 'r-',  'Marker','^', 'linewidth',1.5);
loglog(res1s(:,end,1), res1s(:,10,1), 'g-.', 'Marker','^', 'linewidth',1.5);
loglog(res1s(:,end,2), res1s(:, 4,2), 'b--', 'Marker','o', 'linewidth',1.5);
loglog(res1s(:,end,2), res1s(:, 7,2), 'r-',  'Marker','o', 'linewidth',1.5);
loglog(res1s(:,end,2), res1s(:,10,2), 'g-.', 'Marker','o', 'linewidth',1.5);
loglog(res1s(:,end,3), res1s(:, 4,3), 'b--', 'Marker','s', 'linewidth',1.5);
loglog(res1s(:,end,3), res1s(:, 7,3), 'r-',  'Marker','s', 'linewidth',1.5);
loglog(res1s(:,end,3), res1s(:,10,3), 'g-.', 'Marker','s', 'linewidth',1.5);
loglog(res1s(:,end,4), res1s(:, 4,4), 'b--', 'Marker','d', 'linewidth',1.5);
loglog(res1s(:,end,4), res1s(:, 7,4), 'r-',  'Marker','d', 'linewidth',1.5);
loglog(res1s(:,end,4), res1s(:,10,4), 'g-.', 'Marker','d', 'linewidth',1.5);
loglog(res1s(:,end,5), res1s(:, 4,5), 'b--', 'Marker','v', 'linewidth',1.5);
loglog(res1s(:,end,5), res1s(:, 7,5), 'r-',  'Marker','v', 'linewidth',1.5);
loglog(res1s(:,end,5), res1s(:,10,5), 'g-.', 'Marker','v', 'linewidth',1.5);
legend('OCC-2','POL-2','QFIter-2','OCC-3','POL-3','QFIter-3','OCC-4','POL-4','QFIter-4','OCC-5','POL-5','QFIter-5','OCC-6','POL-6','QFIter-6',...
       'orientation', 'horizontal', 'numcolumns',3, 'location','best');
xlabel('\alpha / °');  ylabel('\epsilon / (°/h)');




